home *** CD-ROM | disk | FTP | other *** search
/ Amiga Format CD 46 / Amiga Format CD46 (1999-10-20)(Future Publishing)(GB)[!][issue 1999-12].iso / -in_the_mag- / reader_requests / scilab / tests / bezout.tst < prev    next >
Text File  |  1999-09-16  |  5KB  |  117 lines
  1. //FICHIER_DE_TEST ADD NAME=bezout.tst,SSI=0
  2. mode(5)
  3. //test
  4. un=poly(1,'s','c');
  5. zer=0*un;
  6. s=poly(0,'s');
  7. [p,q]=bezout(un,s);
  8. if norm(coeff([un s]*q-[p 0]))>10*%eps then pause,end
  9. [p,q]=bezout(s,un);
  10. if norm(coeff([s un]*q-[p 0]))>10*%eps then pause,end
  11. [p,q]=bezout(un,un);
  12. if norm(coeff([un un]*q-[p 0]))>10*%eps then pause,end
  13. //
  14. [p,q]=bezout(zer,s);
  15. if norm(coeff([zer s]*q-[p 0]))>10*%eps then pause,end
  16. [p,q]=bezout(s,zer);
  17. if norm(coeff([s zer]*q-[p 0]))>10*%eps then pause,end
  18. [p,q]=bezout(zer,zer);
  19. if norm(coeff([zer zer]*q-[p 0]))>10*%eps then pause,end
  20. //
  21. [p,q]=bezout(zer,un);
  22. if norm(coeff([zer un]*q-[p 0]))>10*%eps then pause,end
  23. [p,q]=bezout(un,zer);
  24. if norm(coeff([un zer]*q-[p 0]))>10*%eps then pause,end
  25. //
  26.  
  27. //simple
  28. a=poly([1 2 3],'z');
  29. b=poly([4 1],'z');
  30.  
  31. [p q]=bezout(a,b);
  32. dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
  33. if norm(coeff(dt/dt0-1))>10*%eps then pause,end
  34. qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
  35. if norm(coeff([p 0]*qi-[a b]))>100*%eps then pause,end
  36.  
  37. //more difficult
  38.  
  39. b1=poly([4 1+1000*%eps],'z');del=10*norm(coeff(b1-b));
  40. [p,q]=bezout(a,b1);
  41. dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
  42. if norm(coeff(dt/dt0-1))>del then pause,end
  43. qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
  44. if norm(coeff([p 0]*qi-[a b1]))>del then pause,end
  45.  
  46. b1=poly([4 1+.001],'z');del=10*norm(coeff(b1-b));
  47. [p,q]=bezout(a,b1);
  48. dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
  49. if norm(coeff(dt/dt0-1))>100000*%eps then pause,end
  50. qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
  51. if norm(coeff([p 0]*qi-[a b1]))>100000*%eps then pause,end
  52.  
  53.  
  54. b1=poly([4 1+10000*%eps],'z');del=10*norm(coeff(b1-b));
  55. [p,q]=bezout(a,b1);
  56. dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
  57. if norm(coeff(dt/dt0-1))>del then pause,end
  58. qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
  59. if norm(coeff([p 0]*qi-[a b1]))>del then pause,end
  60. //
  61. z=poly(0,'z');
  62. num = 0.99999999999999922+z*(-4.24619123578530730+..
  63.       z*(10.0503215227460350+z*(-14.6836461849931740+..
  64.       z*(13.924822877119892+z*(-5.63165998008533460+..
  65.       z*(-5.63165998008530710+z*(13.9248228771198730+..
  66.       z*(-14.683646184993167+z*(10.0503215227460330+..
  67.       z*(-4.24619123578530910+z*(0.99999999999999989)))))))))));
  68. den = -0.17323463717888873+z*(1.91435457459735380+..
  69.       z*(-9.90126732768255560+z*(31.6286096652930410+..
  70.       z*(-69.3385546880304280+z*(109.586435800377690+..
  71.       z*(-127.516160100808290+z*(109.388684898145950+..
  72.       z*(-67.92645394857864+z*(29.1602681026148110+..
  73.       z*(-7.8212498781094952+z*(0.99999999999999989)))))))))));
  74. //
  75. [p,q]=bezout(num,den);del=1.d-4;
  76. dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
  77. if norm(coeff(dt/dt0-1))>del then pause,end
  78. qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
  79. // JPC 
  80. del=3*del
  81. if norm(coeff([p 0]*qi-[num den]))>del then pause,end
  82. if degree(p)>0 then pause,end
  83. // une autre "solution" telle que l'erreur directe est petite mais l'erreur
  84. // inverse grande
  85. q1=[]
  86. q1(1,1)=poly([0.011533588674 , 3.570224012502 , -28.78817740957 ,...
  87.             102.9479419903, -210.8258579715 , 266.2028963639 ,...
  88.            -207.427018299 , 92.74478640319, -18.5259652457],'z','c');
  89. q1(1,2)=poly([0.000270220664 , -0.002465565223 , 0.010039706635 ,...
  90.            -0.023913827007, 0.036387144032 , -0.036175176058 ,...
  91.             0.022954475171 , -0.008514798968,  0.001417382492],'z','c');
  92. q1(2,1)=poly([0.000639302018 , 20.502472606 , -26.66621972106 ,...
  93.             39.74001534015,  -5.945830824342 , -7.973226036298 ,...
  94.             39.84118622788 , -26.51337424414, 18.5259652457],'z','c');
  95. q1(2,2) =poly( [0.001562930385 , -0.003589174974 , 0.005076237479 ,...
  96.             -0.003483568682,  -0.000135940266 , 0.003651509121 ,...
  97.             -0.005059502869 , 0.003447573440, -0.001417382492],'z','c');
  98. p1 =poly( [0.011422839421 , -0.029264363070 , 0.030070175223 ,...
  99.          -0.012596066108],'z','c');
  100.  
  101. //
  102. //simplification
  103. num =[0.03398330733500143,-0.20716057008572933,0.64660689206696986,...
  104.      -1.97665462134021790,3.38751027286891300,-3.58940006392108120,...
  105.       5.09956159043231680,5.2514918861834694,1.00000000000000020];
  106. den = [0.03398330733500360,-0.20716057008816283,0.64660689204312,...
  107.       -1.97665462079896660,3.38751027286891300,-3.58940006392108350,...
  108.        5.09956159043232040,5.2514918861834712,1];
  109. num=poly(num,'z','c');
  110. den=poly(den,'z','c');
  111. [p,q]=bezout(num,den);del=1.d-8;
  112. dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
  113. if norm(coeff(dt/dt0-1))>del then pause,end
  114. qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
  115. if norm(coeff([p 0]*qi-[num den]))>del then pause,end
  116. if degree(p)<8 then pause,end
  117.